Checking past and future climates at the location of the populations

Published

April 30, 2024

1 Introduction

Goal: compare the values of the climatic variables at the location of the populations under current and future climates. The climatic data can be generated in two ways:

  • point estimates at the location of the population, generated using scale-free downscaling, resulting in higher precision.

  • extracted values from rasters. The pixel values in the raster file correspond to the mean values of the climatic variables across the period considered.

Comment: in all the plots below, the points corresponds to the 34 studied populations of maritime pine and the colors correspond to the main gene pool they belong to.

Climatic data comes from the Climate Downscaling Tool (ClimateDT) described in Marchi et al. (2024).

Set of the climatic variables of interest:

Code
# Selected climatic variables
# ===========================
clim_var <- readRDS(here("data/ClimaticData/NamesSelectedVariables.rds"))

extract_climatedt_metadata(var_clim = clim_var) %>% 
  dplyr::select(label,description,unit) %>% 
  set_colnames(str_to_title(colnames(.))) %>% 
  kable_mydf(font_size = 12)
Label Description Unit
bio1 Mean annual temperature Celsius degrees (°C)
bio3 Isothermality (bio2/bio7) (×100) index
bio4 Temperature seasonality (standard deviation ×100) Celsius degrees (°C)
bio12 Annual precipitation Millimeters (mm)
bio15 Precipitation seasonality (coefficient of variation) index
SHM Summer heat moisture index °C / mm

2 Comparing past reference periods

To account for the climates under which the populations have evolved, we could use the reference period 1901-1950 or 1961-1990.

Until now, we have used the reference period 1901-1950 to capture as much as possible the climate under which the populations have evolved, i.e. before the human-induced climate change. Note that the human induced climate change may be minor before 1990 so we may use the pre-industrialization period 1961-1990 period, which may not be very different from the 1901-1950 period.

Another reason we used the reference period 1901-1950 is that our reference period has to be before 1983 because we have to use the gene-climate relationships estimated under this reference period to then predict the disruption of the gene-environment relationships under the climates between 1983 and 2014 at the location of the National Forest Inventory plots (one of the validation steps of the genomic offset method we do in the paper).

Note that the historic climatic layers for the 1901-1950 period may not be very accurate according to Maurizio Marchi, that’s why he suggested to use the 1961-1990 period.

Below, we compare the mean values of the climatic variables of interest for the reference period 1901-1950 and 1961-1990. We use point estimates climatic data adjusted or not for elevation.

Code
# Load population past climate (point estimates) in a list with data adjusted or not for elevation
# The rds file corresponds to a list with different reference periods, here 1901-1950 and 1961-1990
list_clim_past <- list(ADJ=list(name="adjusted for elevation",
                                code="ADJ",
                                data=readRDS(here("data/ClimaticData/MaritimePinePops/ClimatePopulationLocationPointEstimates_ReferencePeriods_ADJ.rds"))),
                       noADJ=list(name="NOT adjusted for elevation",
                                  code="noADJ",
                                  data=readRDS(here("data/ClimaticData/MaritimePinePops/ClimatePopulationLocationPointEstimates_ReferencePeriods_noADJ.rds"))))
 
# Load information about the main gene pool of each population (for the figures)
gps <- readRDS(here("data/GenomicData/MainGenePoolPopulations.rds")) %>% arrange(pop)

# Population coordinates as spatial points (for the extraction of future climates)
pop_coord <-  list_clim_past[[1]]$data$ref_1901_1950$ref_means %>% 
  dplyr::select(longitude,latitude) %>% 
  SpatialPoints(proj4string=CRS("+proj=longlat +datum=WGS84 +no_defs"))
Code
# Function to generate the scatter plots
# ======================================
make_ggscatter_plots <- function(df,x,gps,xlab,ylab){
  
var_name <- extract_climatedt_metadata(var_clim = x) %>% 
  mutate(var_legend= paste0(description, " (", label,"; ",unit_symbol,")")) %>% 
  pull(var_legend)
    
df %>% 
    left_join(gps, by="pop") %>% 
    mutate(main_gp_pop=factor(main_gp_pop,levels=c("Northern Africa",
                                                   "Corsica",
                                                   "Central Spain",
                                                   "French Atlantic",
                                                   "Iberian Atlantic",
                                                   "South-eastern Spain"))) %>% 
    
  ggplot(aes(x=x,y=y)) + 
  geom_point(aes(color=main_gp_pop),size=3) +
  scale_colour_manual(name="Main gene pool",
                      values = c("orangered3","gold2","darkorchid3","navyblue","turquoise2","green3")) +
  geom_abline(intercept = 0, slope = 1, color="gray60") +
  xlab(xlab) +
  ylab(ylab) +
  ggtitle(var_name) + 
  theme_bw() +
  theme(legend.position = "none",
        legend.title = element_text(size=16, color=bg_color),
        plot.title = element_text(size=16),
        legend.text = element_text(size = 16, color=bg_color),
        axis.text = element_text(size = 14, color=bg_color),
        axis.title = element_text(size=16, color=bg_color))
    
  
  }
Code
lapply(list_clim_past, function(clim_past){
  
  
# Generate the scatter plots
scatter_plots <- lapply(clim_var, function(x){

  lapply(clim_past$data, function(y){
    
    y$ref_means %>% 
      dplyr::select(pop,all_of(x))
  }) %>% 
    bind_rows(.id="ref_period") %>% 
    pivot_wider(names_from = ref_period,values_from = x) %>% 
    set_colnames(c("pop","x","y"))  %>% 
    make_ggscatter_plots(x=x,
                       gps=gps,
                       xlab=paste0("Reference period ",clim_past$data[[1]]$range[[1]],"-",clim_past$data[[1]]$range[[2]]),
                       ylab=paste0("Reference period ",clim_past$data[[2]]$range[[1]],"-",clim_past$data[[2]]$range[[2]]))
})


p <- plot_grid(plotlist=scatter_plots)

# legend 
legend <- get_legend(scatter_plots[[1]] + theme(legend.position = c(0.5,0.5)) )

# merge plots + legend
p <- plot_grid(p, legend, ncol = 2,rel_widths = c(1, 0.2))


# make title
title <- ggdraw() + 
  draw_label(paste0("Climatic data ",clim_past$name),
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))
  
  
})

Interpretation:

  • Similar values between adjusted or not adjusted climatic data.

  • The mean climatic values across the reference periods 1901-1950 and 1961-1990 are mostly similar for AHM and SHM.

  • Increase in mean annual temperature (bio1) for all populations across the 1961-1990 period compared to the 1901-1950 period.

  • Populations from the French Atlantic gene pool experienced a decrease in isothermality and Hargreaves reference evaporation and an increase in mean temperature of the driest quarter.

  • Populations from the Iberian gene pool experienced an increase in isothermality, mean annual temperature and Hargreaves reference evapotranspiration.

3 Point estimates vs raster extraction

Code
# Function to extract the climatic variables from rasters
# =======================================================

extract_clim_var <- function(x,gcm,period,ssp){
  
if(period=="1901-1950"){
  path <- paste0("data/ClimaticData/ClimateDTRasters/1km_1901-1950_Extent-JulietteA/",x,".tif")} else {
  path <- paste0("data/ClimaticData/ClimateDTRasters/1km_",gcm,"_",period,"_",ssp,"_Extent-JulietteA/",x,".tif")
  }
  
# For the variable Eref, there are some NAs for some populations
# So we start with a simple extraction (extraction of the cell value if the point falls within the cell)
# and if there are some NAs, we perform a bilinear extraction
  # bilinear interpolation = the returned values are interpolated from the values of the four nearest raster cells.
# And we replace the NAs by the extracted values obtained with the bilinear method
# if there are still NAs, we extract cell values using a buffer of 2km
  # If the distance between the sampling point and the center of a cell is less than or equal to the buffer, the cell is included. 
# We take the mean (fun=mean) of the cell values obtained with the extraction done with the buffer method
# And we replace the NAs by the values obtained with the buffer method (with a buffer of 2km)
# If there are still NAs, we do the same with a buffer of 3km
  
ext_simple <- raster::raster(here(path)) %>% raster::extract(pop_coord)
i <- is.na(ext_simple)

if(length(which(i==TRUE)) !=0){
  
ext_bilinear <- raster::raster(here(path)) %>% raster::extract(pop_coord,method="bilinear")
ext_simple[i] <- ext_bilinear[i]
i <- is.na(ext_simple)

if(length(which(i==TRUE)) !=0){
  
ext_buff <- raster::raster(here(path)) %>% raster::extract(pop_coord,buffer=2000,fun=mean)
ext_simple[i] <- ext_buff[i]
i <- is.na(ext_simple)

if(length(which(i==TRUE)) !=0){
  
ext_buff <- raster::raster(here(path)) %>% raster::extract(pop_coord,buffer=3000,fun=mean)
ext_simple[i] <- ext_buff[i]
i <- is.na(ext_simple)


}

}
  
}

return(ext_simple)

}

3.1 Past climate (1901-1950)

We compare the point estimates climatic values (from csv) with extracted values from raster files for the period 1901-1950.

Code
lapply(list_clim_past, function(clim_past){
  
  
scatter_plots <- lapply(clim_var, function(x) {
  
    
clim_past$data$ref_1901_1950$ref_means[,c("pop",x)] %>% 
    mutate(ext_var =extract_clim_var(x=x,period="1901-1950")) %>% 
    set_colnames(c("pop","x","y")) %>% 
    make_ggscatter_plots(x=x,
                       gps=gps,
                       xlab="Point estimates",
                       ylab="Extracted values from raster")
  }) 

p <- plot_grid(plotlist=scatter_plots)

# legend 
legend <- get_legend(scatter_plots[[1]] + theme(legend.position = c(0.5,0.5)) )

# merge plots + legend
p <- plot_grid(p, legend, ncol = 2,rel_widths = c(1, 0.2))


# make title
title <- ggdraw() + 
  draw_label(paste0("Point estimates ",clim_past$name),
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))
  
  
})

3.2 Future climate (2041-2070 - GFDL-ESM4 - SSP 3.7-0)

Code
# Point estimate future climatic values
list_fut_clim <-  list(ADJ=list(name="adjusted for elevation",
                                code="ADJ",
                                data=readRDS(here("data/ClimaticData/MaritimePinePops/ClimatePopulationLocationsPointEstimates_GFDL-ESM4_SSP370_2041-2070.rds"))[[1]]),
                       noADJ=list(name="NOT adjusted for elevation",
                                  code="noADJ",
                                  data=readRDS(here("data/ClimaticData/MaritimePinePops/ClimatePopulationLocationsPointEstimates_GFDL-ESM4_SSP370_2041-2070.rds"))[[2]]))                    
                      

lapply(list_fut_clim, function(fut_clim){
  
  
scatter_plots <- lapply(clim_var, function(x) {
  
    
fut_clim$data[,c("pop",x)] %>% 
    mutate(ext_var =extract_clim_var(x=x,period="2041-2070",gcm="GFDL-ESM4",ssp="ssp370")) %>% 
    set_colnames(c("pop","x","y")) %>% 
    make_ggscatter_plots(x=x,
                       gps=gps,
                       xlab="Point estimates",
                       ylab="Extracted values from raster")
  }) 

p <- plot_grid(plotlist=scatter_plots)

# legend 
legend <- get_legend(scatter_plots[[1]] + theme(legend.position = c(0.5,0.5)) )

# merge plots + legend
p <- plot_grid(p, legend, ncol = 2,rel_widths = c(1, 0.2))


# make title
title <- ggdraw() + 
  draw_label(paste0("Point estimates ",fut_clim$name),
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))
  
  
})

Conclusion: Point estimate values and climatic values extracted from rasters are very similar, thus suggesting that we can use them interchangeably. Point estimate values adjusted or not for elevation are also similar, so using one or the other should not influence the results of the following analyses.

4 Comparing past and future climates

We compare the past and future climates using either:

  • the raw difference between the mean future climate and mean past climate (i.e. mean future climate - mean past climate)

  • the relative climatic distance: relative distance = ( future climatic value - past climatic value ) / past climatic value.

4.1 Raw differences between past and future climates

We calculate the raw difference between the mean future climate and mean past climate (i.e. mean future climate - mean past climate) using:

  • climatic values extracted from rasters for the future climates of the five global climate models (GCMs) for the period 2041-2070 and SSP3-7.0.

  • past climatic values either from point estimates or extracted from rasters (reference period 1901-1950).

Code
# Point estimate future climatic values (2041-2070 - GFDL-ESM4 - SSP 3.7-0)
fut_clim <- readRDS(here("data/ClimaticData/MaritimePinePops/ClimatePopulationLocationsPointEstimates_GFDL-ESM4_SSP370_2041-2070.rds"))

# GCM names
gcm_names <- list.files(here("data/ClimaticData/ClimateDTRasters/")) %>%  str_sub(5,-35) %>% setdiff("")

# Extract population future climate for each GCM and merge with past climatic data
fut_rast <- lapply(list_clim_past, function(clim_past){
  
lapply(clim_past$data, function(y){

  
df_ref_means <-  y$ref_means %>% 
  dplyr::select(pop,all_of(clim_var)) %>% 
  pivot_longer(cols=clim_var,names_to="variable",values_to="mean_ref")
  
lapply(clim_var, function(x){
  
lapply(gcm_names, function(gcm) extract_clim_var(x=x,gcm=gcm,period="2041-2070",ssp="ssp370")) %>% 
  setNames(gcm_names) %>% 
  as_tibble() %>% 
  mutate(pop=unique(clim_past$data[[1]]$ref_means$pop))

}) %>% setNames(clim_var) %>% 
  bind_rows(.id = "variable") %>% 
  pivot_longer(names_to = "gcm",values_to="mean_fut",cols=all_of(gcm_names)) %>% 
  mutate(fut_method= "extracted from rasters") %>% 
  mutate(ref_period=paste0(y$range[[1]],"-",y$range[[2]])) %>% 
  left_join(df_ref_means, by=c("pop","variable")) %>% 
  left_join(gps, by="pop") 
    

}) %>% 
    bind_rows() %>% 
    mutate(adj=clim_past$code)
  
}) %>%  bind_rows()
Warning: Using an external vector in selections was deprecated in tidyselect 1.1.0.
ℹ Please use `all_of()` or `any_of()` instead.
  # Was:
  data %>% select(clim_var)

  # Now:
  data %>% select(all_of(clim_var))

See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
Code
violin_plots <- lapply(unique(fut_rast$adj), function(y){

# We generate a table of 1020 rows: 34 populations * 5 GCMs * 6 climatic variables
clim_df <- fut_rast %>% 
  filter(ref_period =="1901-1950") %>% 
  filter(adj==y)


clim_df <- lapply(clim_var, function(x) {
  
    
list_clim_past$ADJ$data$ref_1901_1950$ref_means[,"pop"] %>% 
    mutate(mean_ref_rast =extract_clim_var(x=x,period="1901-1950"),
           variable=x) 
  }) %>% 
  bind_rows() %>% 
  right_join(clim_df,by = c("pop", "variable"))



clim_df %>% 
  pivot_longer(cols=c("mean_ref_rast","mean_ref"),names_to="mean_ref_meth",values_to = "mean_ref") %>% 
  group_split(mean_ref_meth) %>% 
  lapply(function(meth){


violin_plots <- lapply(clim_var, function(x) {
  
var_name <- extract_climatedt_metadata(var_clim = x) %>% 
  mutate(var_legend= paste0(description, " (", label,"; ",unit_symbol,")")) %>% 
  pull(var_legend)
  
meth %>% 
    dplyr::filter(variable %in% x) %>% 
    mutate(diff=mean_fut-mean_ref) %>% 
    mutate(main_gp_pop=factor(main_gp_pop,levels=c("Northern Africa",
                                                 "Corsica",
                                                 "Central Spain",
                                                 "French Atlantic",
                                                 "Iberian Atlantic",
                                                 "South-eastern Spain"))) %>% 
  ggplot(aes(x=gcm, y=diff)) + 
  geom_hline(yintercept=0, color=bg_color) +
  geom_violin(alpha=0.2, fill="gray76", color=bg_color) + 
  geom_jitter(shape=16,aes(colour=main_gp_pop),position=position_jitter(0.2),size=2) +
  scale_colour_manual(values = c("orangered3",
                                                             "gold2",
                                                             "darkorchid3",
                                                              "navyblue",
                                                             "turquoise2",
                                                             "green3")) +
  xlab("") + 
  ylab("mean future climate - mean past climate") +
  labs(colour="Main gene pool") +
  theme_bw() +
  scale_x_discrete(guide = guide_axis(angle = 90)) +
  ggtitle(var_name) +
  theme(legend.position="none",
        legend.title = element_text(size=14, color=bg_color),
        legend.text = element_text(size = 14, color=bg_color),
        axis.text = element_text(size = 14, color=bg_color), 
        axis.title = element_text(size=16, color=bg_color),
        ) + 
  guides(colour = guide_legend(override.aes = list(size=4))) 
})

p <- plot_grid(plotlist=violin_plots)

# legend 
legend <- get_legend(violin_plots[[6]] + theme(legend.position = c(0.5,0.5)) )

# merge plots + legend
p <- plot_grid(p, legend, ncol = 2,rel_widths = c(1, 0.2))


# make title
if(unique(meth$mean_ref_meth)=="mean_ref"){
  name_title <- paste0("Past climatic values from point estimates ", list_clim_past[[y]]$name, " - reference period 1901-1950")
} else {
  name_title <- paste0("Past climatic values extracted from raster files - reference period 1901-1950")
}

title <- ggdraw() + 
  draw_label(name_title,
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))

})

})


violin_plots[[1]]

Code
violin_plots[[2]][[1]]

Predictions consistent across GCMs.

4.2 Relative climatic distance

To compare past and future climatic values, we calculate, for each GCM and each climatic variable, the relative distance between past and future climatic values, i.e relative distance = ( future climatic value - past climatic value ) / past climatic value.

4.2.1 Future climatic data from rasters

On the graph below, the relative distances are on the y-axis and are calculated using point estimates for the past climates (i.e. climates under the reference periods 1901-1950 or 1961-1990) and raster data for the future climates (mean predicted climatic values across the period 2041-2070 and the SSP3-7.0. The past climates correspond to:

  • either the reference period 1901-1950 or 1961-1990.

  • either the adjusted or not adjusted data for elevation.

Code
violin_plots <- lapply(unique(fut_rast$adj), function(y){
  
lapply(unique(fut_rast$ref_period), function(z){

violin_plots <- lapply(gcm_names, function(x) {
  
fut_rast %>% 
  dplyr::filter(gcm %in% x) %>% 
  dplyr::filter(ref_period %in% z) %>% 
  dplyr::filter(adj %in% y) %>% 
  mutate(rel_diff=(mean_fut-mean_ref)/mean_ref) %>% 
  mutate(main_gp_pop=factor(main_gp_pop,levels=c("Northern Africa",
                                                 "Corsica",
                                                 "Central Spain",
                                                 "French Atlantic",
                                                 "Iberian Atlantic",
                                                 "South-eastern Spain"))) %>% 
  ggplot(aes(x=variable, y=rel_diff)) + 
  geom_hline(yintercept=0, color=bg_color) +
  geom_violin(alpha=0.2, fill="gray76", color=bg_color) + 
  geom_jitter(shape=16,aes(colour=main_gp_pop),position=position_jitter(0.2),size=1.7) +
  scale_colour_manual(values = c("orangered3","gold2","darkorchid3","navyblue","turquoise2","green3")) +
  xlab("") + 
  ylab("Relative climatic distance") +
  labs(colour="Main gene pool") +
  theme_bw() +
  ggtitle(x) +
  theme(legend.position="none",
        legend.title = element_text(size=14, color=bg_color),
        legend.text = element_text(size = 14, color=bg_color),
        axis.text = element_text(size = 14, color=bg_color),
        axis.title = element_text(size=16, color=bg_color)) + 
  guides(colour = guide_legend(override.aes = list(size=4))) 
})


# legend 
violin_plots$legend <- get_legend(violin_plots[[1]] + theme(legend.position = c(0.5,0.5)) )


# make title
title <- ggdraw() + 
  draw_label(paste0("Future climate from rasters vs reference period ",z," (point estimates ",list_clim_past[[y]]$name,")"),
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots
p <- plot_grid(plotlist=violin_plots)

# ========================================
# Figure for the Supplementary Information

# We save plots for the reference period 1901-1950 and with point estimates not adjusted for elevation
if (y=="noADJ" & z=="1901-1950"){

 ggsave(p,device="pdf",filename=here("figs/ExploratoryAnalyses/RelativeClimaticDistances_SelectedVariables.pdf"),
        width=12,height=10)
    
}
# ========================================

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))
          
  
})

})
Warning in get_plot_component(plot, "guide-box"): Multiple components found; returning the first one. To return all, use `return_all = TRUE`.
Warning in get_plot_component(plot, "guide-box"): Multiple components found; returning the first one. To return all, use `return_all = TRUE`.
Warning in get_plot_component(plot, "guide-box"): Multiple components found; returning the first one. To return all, use `return_all = TRUE`.
Warning in get_plot_component(plot, "guide-box"): Multiple components found; returning the first one. To return all, use `return_all = TRUE`.
Code
violin_plots

We also create a figure for the Supplementary Information in which we highlight in pink the populations ALT and ARM, i.e. the two populations considered the most at risk of maladaptation with GDM and GF based on the candidate SNPs (the two methods for which GO predictions were the best predictors of mortality rates in the NFI populations and the common gardens).

Code
y="noADJ" # point estimates not adjusted for elevation
z="1901-1950" # we use the reference period 1901-1950 

violin_plots <- lapply(gcm_names, function(x) {

df <- fut_rast %>% 
  dplyr::filter(gcm %in% x) %>% 
  dplyr::filter(ref_period %in% z) %>% 
  dplyr::filter(adj %in% y) %>% 
  mutate(rel_diff=(mean_fut-mean_ref)/mean_ref) %>% 
  mutate(pop_to_color=ifelse(pop %in% c("ALT","ARM"), pop, "Other populations"))

  ggplot(df,aes(x=variable, y=rel_diff)) + 
  geom_hline(yintercept=0, color=bg_color) +
  geom_violin(alpha=0.2, fill="gray76", color=bg_color) + 
  geom_jitter(data=df %>% dplyr::filter(pop_to_color %in% "Other populations"),shape=16,colour="lightskyblue1",position=position_jitter(0.2),size=1.7) +
  geom_jitter(data=df %>% dplyr::filter(pop_to_color %in% c("ALT","ARM")),shape=16,colour="deeppink",position=position_jitter(0.2),size=2) +
  xlab("") + 
  ylab("Relative climatic distance") +
  labs(colour="Main gene pool") +
  theme_bw() +
  ggtitle(x) +
  theme(legend.position="none",
        legend.title = element_text(size=14, color=bg_color),
        legend.text = element_text(size = 14, color=bg_color),
        axis.text = element_text(size = 14, color=bg_color),
        axis.title = element_text(size=16, color=bg_color)) + 
  guides(colour = guide_legend(override.aes = list(size=4))) 
})


# merge plots
p <- plot_grid(plotlist=violin_plots)

# We save the plot for the Supplementary Information:

 ggsave(p,
        device="pdf",
        filename=here("figs/ExploratoryAnalyses/RelativeClimaticDistances_SelectedVariables_2PopsHighlighted.pdf"),
        width=12,
        height=10)
    
p

Interpretation: For all populations, AHM, SHM, bio1 and bio9 and Eref are predicted to increase under future climates. Some exceptions for AHM for some populations from the Iberian and French Atlantic regions for the MRI-ESM2-0 GCM. Different predictions across GCMs for the OLB populations. The different predictions across GCMs show the importance of calculating the genomic offset under the different GCMs.

4.2.2 Future climatic data from point estimates

In the graphs below, the future and past climates correspond to point estimates adjusted or not for elevation.

Code
violin_plots <- lapply(unique(fut_rast$adj), function(y){
  
violin_plots <- lapply(unique(fut_rast$ref_period), function(z){
  

df <- fut_rast %>% 
  dplyr::filter(ref_period %in% z) %>% 
  dplyr::filter(adj %in% y) %>% 
  dplyr::filter(gcm %in% "GFDL-ESM4") %>% 
  dplyr::select(-mean_fut, -fut_method,-adj)
  
fut_clim[[y]] %>%
  dplyr::select(pop,all_of(clim_var)) %>% 
  pivot_longer(cols=all_of(clim_var), names_to="variable",values_to="mean_fut") %>% 
  arrange(variable,pop) %>% 
  mutate(gcm="GFDL-ESM4") %>% 
  inner_join(df, by=c("pop","variable","gcm")) %>% 
  dplyr::mutate(rel_diff=(mean_fut-mean_ref)/mean_ref) %>% 
  
  mutate(main_gp_pop=factor(main_gp_pop,levels=c("Northern Africa",
                                                 "Corsica",
                                                 "Central Spain",
                                                 "French Atlantic",
                                                 "Iberian Atlantic",
                                                 "South-eastern Spain"))) %>% 
  ggplot(aes(x=variable, y=rel_diff)) + 
  geom_hline(yintercept=0, color=bg_color) +
  geom_violin(alpha=0.2, fill="gray76", color=bg_color) + 
  geom_jitter(shape=16,aes(colour=main_gp_pop),position=position_jitter(0.2),size=1.7) +
  scale_colour_manual(values = c("orangered3",
                                                             "gold2",
                                                             "darkorchid3",
                                                              "navyblue",
                                                             "turquoise2",
                                                             "green3")) +
  xlab("") + 
  ylab("Relative climatic distance") +
  labs(colour="Main gene pool") +
  theme_bw() +
  ggtitle(paste0("Climatic data ",list_clim_past[[y]]$name," - reference period ",z)) +
  theme(legend.position="none",
        legend.title = element_text(size=14, color=bg_color),
        legend.text = element_text(size = 14, color=bg_color),
        axis.text = element_text(size = 14, color=bg_color),
        axis.title = element_text(size=16, color=bg_color)) + 
  guides(colour = guide_legend(override.aes = list(size=4))) 
  
})
  
plot_grid(plotlist=violin_plots)
  
}) 


# make title
title <- ggdraw() + 
  draw_label("Future and past climate from point estimates - GCM GFDL-ESM4",
             fontface = 'bold',
             x = 0,
             hjust = 0) +
  theme(plot.margin = margin(0, 0, 0, 7))

# merge plots
p <- plot_grid(plotlist=violin_plots,nrow=2)

# merge plots + title
plot_grid(title, p, ncol = 1,rel_heights = c(0.1, 1))

5 Session information

Code
devtools::session_info()
─ Session info ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
 setting  value
 version  R version 4.4.0 (2024-04-24)
 os       Ubuntu 22.04.4 LTS
 system   x86_64, linux-gnu
 ui       X11
 language (EN)
 collate  fr_FR.UTF-8
 ctype    fr_FR.UTF-8
 tz       Europe/Paris
 date     2024-05-01
 pandoc   3.1.11 @ /usr/lib/rstudio/resources/app/bin/quarto/bin/tools/x86_64/ (via rmarkdown)

─ Packages ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
 package     * version date (UTC) lib source
 cachem        1.0.8   2023-05-01 [1] CRAN (R 4.4.0)
 cellranger    1.1.0   2016-07-27 [1] CRAN (R 4.4.0)
 cli           3.6.2   2023-12-11 [1] CRAN (R 4.4.0)
 codetools     0.2-20  2024-03-31 [4] CRAN (R 4.4.0)
 colorspace    2.1-0   2023-01-23 [1] CRAN (R 4.4.0)
 cowplot     * 1.1.3   2024-01-22 [1] CRAN (R 4.4.0)
 devtools      2.4.5   2022-10-11 [1] CRAN (R 4.4.0)
 digest        0.6.35  2024-03-11 [1] CRAN (R 4.4.0)
 dplyr       * 1.1.4   2023-11-17 [1] CRAN (R 4.4.0)
 ellipsis      0.3.2   2021-04-29 [1] CRAN (R 4.4.0)
 evaluate      0.23    2023-11-01 [1] CRAN (R 4.4.0)
 fansi         1.0.6   2023-12-08 [1] CRAN (R 4.4.0)
 farver        2.1.1   2022-07-06 [1] CRAN (R 4.4.0)
 fastmap       1.1.1   2023-02-24 [1] CRAN (R 4.4.0)
 forcats     * 1.0.0   2023-01-29 [1] CRAN (R 4.4.0)
 fs            1.6.4   2024-04-25 [1] CRAN (R 4.4.0)
 generics      0.1.3   2022-07-05 [1] CRAN (R 4.4.0)
 ggplot2     * 3.5.1   2024-04-23 [1] CRAN (R 4.4.0)
 glue          1.7.0   2024-01-09 [1] CRAN (R 4.4.0)
 gtable        0.3.5   2024-04-22 [1] CRAN (R 4.4.0)
 here        * 1.0.1   2020-12-13 [1] CRAN (R 4.4.0)
 highr         0.10    2022-12-22 [1] CRAN (R 4.4.0)
 hms           1.1.3   2023-03-21 [1] CRAN (R 4.4.0)
 htmltools     0.5.8.1 2024-04-04 [1] CRAN (R 4.4.0)
 htmlwidgets   1.6.4   2023-12-06 [1] CRAN (R 4.4.0)
 httpuv        1.6.15  2024-03-26 [1] CRAN (R 4.4.0)
 jsonlite      1.8.8   2023-12-04 [1] CRAN (R 4.4.0)
 kableExtra  * 1.4.0   2024-01-24 [1] CRAN (R 4.4.0)
 knitr       * 1.46    2024-04-06 [1] CRAN (R 4.4.0)
 labeling      0.4.3   2023-08-29 [1] CRAN (R 4.4.0)
 later         1.3.2   2023-12-06 [1] CRAN (R 4.4.0)
 lattice       0.22-6  2024-03-20 [4] CRAN (R 4.4.0)
 lifecycle     1.0.4   2023-11-07 [1] CRAN (R 4.4.0)
 lubridate   * 1.9.3   2023-09-27 [1] CRAN (R 4.4.0)
 magrittr    * 2.0.3   2022-03-30 [1] CRAN (R 4.4.0)
 memoise       2.0.1   2021-11-26 [1] CRAN (R 4.4.0)
 mime          0.12    2021-09-28 [1] CRAN (R 4.4.0)
 miniUI        0.1.1.1 2018-05-18 [1] CRAN (R 4.4.0)
 munsell       0.5.1   2024-04-01 [1] CRAN (R 4.4.0)
 pillar        1.9.0   2023-03-22 [1] CRAN (R 4.4.0)
 pkgbuild      1.4.4   2024-03-17 [1] CRAN (R 4.4.0)
 pkgconfig     2.0.3   2019-09-22 [1] CRAN (R 4.4.0)
 pkgload       1.3.4   2024-01-16 [1] CRAN (R 4.4.0)
 plyr          1.8.9   2023-10-02 [1] CRAN (R 4.4.0)
 profvis       0.3.8   2023-05-02 [1] CRAN (R 4.4.0)
 promises      1.3.0   2024-04-05 [1] CRAN (R 4.4.0)
 purrr       * 1.0.2   2023-08-10 [1] CRAN (R 4.4.0)
 R6            2.5.1   2021-08-19 [1] CRAN (R 4.4.0)
 ragg          1.3.0   2024-03-13 [1] CRAN (R 4.4.0)
 raster      * 3.6-26  2023-10-14 [1] CRAN (R 4.4.0)
 Rcpp          1.0.12  2024-01-09 [1] CRAN (R 4.4.0)
 readr       * 2.1.5   2024-01-10 [1] CRAN (R 4.4.0)
 readxl      * 1.4.3   2023-07-06 [1] CRAN (R 4.4.0)
 remotes       2.5.0   2024-03-17 [1] CRAN (R 4.4.0)
 reshape2    * 1.4.4   2020-04-09 [1] CRAN (R 4.4.0)
 rlang         1.1.3   2024-01-10 [1] CRAN (R 4.4.0)
 rmarkdown     2.26    2024-03-05 [1] CRAN (R 4.4.0)
 rprojroot     2.0.4   2023-11-05 [1] CRAN (R 4.4.0)
 rstudioapi    0.16.0  2024-03-24 [1] CRAN (R 4.4.0)
 scales        1.3.0   2023-11-28 [1] CRAN (R 4.4.0)
 sessioninfo   1.2.2   2021-12-06 [1] CRAN (R 4.4.0)
 shiny         1.8.1.1 2024-04-02 [1] CRAN (R 4.4.0)
 sp          * 2.1-4   2024-04-30 [1] CRAN (R 4.4.0)
 stringi       1.8.3   2023-12-11 [1] CRAN (R 4.4.0)
 stringr     * 1.5.1   2023-11-14 [1] CRAN (R 4.4.0)
 svglite       2.1.3   2023-12-08 [1] CRAN (R 4.4.0)
 systemfonts   1.0.6   2024-03-07 [1] CRAN (R 4.4.0)
 terra         1.7-71  2024-01-31 [1] CRAN (R 4.4.0)
 textshaping   0.3.7   2023-10-09 [1] CRAN (R 4.4.0)
 tibble      * 3.2.1   2023-03-20 [1] CRAN (R 4.4.0)
 tidyr       * 1.3.1   2024-01-24 [1] CRAN (R 4.4.0)
 tidyselect    1.2.1   2024-03-11 [1] CRAN (R 4.4.0)
 tidyverse   * 2.0.0   2023-02-22 [1] CRAN (R 4.4.0)
 timechange    0.3.0   2024-01-18 [1] CRAN (R 4.4.0)
 tzdb          0.4.0   2023-05-12 [1] CRAN (R 4.4.0)
 urlchecker    1.0.1   2021-11-30 [1] CRAN (R 4.4.0)
 usethis       2.2.3   2024-02-19 [1] CRAN (R 4.4.0)
 utf8          1.2.4   2023-10-22 [1] CRAN (R 4.4.0)
 vctrs         0.6.5   2023-12-01 [1] CRAN (R 4.4.0)
 viridisLite   0.4.2   2023-05-02 [1] CRAN (R 4.4.0)
 withr         3.0.0   2024-01-16 [1] CRAN (R 4.4.0)
 xfun          0.43    2024-03-25 [1] CRAN (R 4.4.0)
 xml2          1.3.6   2023-12-04 [1] CRAN (R 4.4.0)
 xtable      * 1.8-4   2019-04-21 [1] CRAN (R 4.4.0)
 yaml          2.3.8   2023-12-11 [1] CRAN (R 4.4.0)

 [1] /home/juliette/R/x86_64-pc-linux-gnu-library/4.4
 [2] /usr/local/lib/R/site-library
 [3] /usr/lib/R/site-library
 [4] /usr/lib/R/library

──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────

References

Marchi, Maurizio, Gabriele Bucci, Paolo Iovieno, and Duncan Ray. 2024. “ClimateDT: A Global Scale-Free Dynamic Downscaling Portal for Historic and Future Climate Data.” Environments 11 (4): 82.